Introduction à l'apprentissage par renforcement¶
TP 1 - les manchots multi-bras¶
1/4 de la note finale est liée à la mise en forme :
- pensez à nettoyer les outputs inutiles (installation, messages de débuggage, ...)
- soignez vos figures : les axes sont-ils faciles à comprendre ? L'échelle est adaptée ?
- commentez vos résultats : vous attendiez-vous à les avoir ? Est-ce étonnant ? Faites le lien avec la théorie.
Ce TP reprend l'exemple d'un médecin et de ses vaccins. Vous allez comparer plusieurs stratégies et trouver celle optimale. Un TP se fait seul ou en binôme. Aucun groupe de plus de 2 personnes.
Vous allez rendre le TP depuis un lien GitHub avec ce notebook mais une version du rapport exportée en PDF & HTML.
%%capture
!pip install matplotlib tqdm numpy ipympl opencv-python
!jupyter labextension install @jupyter-widgets/jupyterlab-manager
!jupyter labextension install jupyter-matplotlib
%load_ext autoreload
%autoreload 2
%matplotlib inline
import typing as t
import math
import torch
import numpy as np
from tqdm.auto import trange, tqdm
import matplotlib.pyplot as plt
from matplotlib.animation import FuncAnimation
import matplotlib.animation as animation
from matplotlib.backends.backend_agg import FigureCanvasAgg as FigureCanvas
import cv2
from IPython.display import display, clear_output, HTML
import random
torch.random.manual_seed(0)
K = 5 # num arms
Présentation du problème¶
class ArmBernoulli:
def __init__(self, p: float):
"""
Vaccine treatment following a Bernoulli law (mean is p and variance is p(1-p)
Args:
p (float): mean parameter
>>> torch.random.manual_seed(random_state)
>>> arm = ArmBernoulli(0.5)
>>> arm.sample(5)
tensor([ True, False, True, True, True])
"""
self.immunity_rate = p
def sample(self, n: int = 1):
return torch.rand(n) < self.immunity_rate
def __repr__(self):
return f'<ArmBernoulli p={self.immunity_rate}'
def generate_arms(num_arms: int):
means = torch.rand(K)
MAB = [ArmBernoulli(m) for m in means]
assert MAB[0].immunity_rate == means[0]
assert (MAB[0].sample(10) <= 1).all() and (MAB[0].sample(10) >= 0).all()
return MAB
MAB = generate_arms(K)
Ce TP reprend l'exemple du médecin présenté en cours.
Q1. Créez une fonction pour trouver $\mu^*$ à partir d'un MAB. Comment est définie la récompense $R_k$ ? Que représente concrètement le regret dans le contexte de ce TP ?
def find_mu_star_index(MAB: t.List[ArmBernoulli]):
return torch.argmax(torch.tensor([arm.immunity_rate for arm in MAB]))
mu_star_index = find_mu_star_index(MAB)
print(f"L'index de mu* est: {mu_star_index}")
def find_mu_star_value(MAB: t.List[ArmBernoulli]):
return max(arm.immunity_rate for arm in MAB)
mu_star_value = find_mu_star_value(MAB)
print(f"La valeur de mu* est: {mu_star_value}")
L'index de mu* est: 1 La valeur de mu* est: 0.7682217955589294
Le regret est une mesure qui représente la différence entre ce que le médecin aurait pu obtenir s'il avait toujours choisi le meilleur vaccin et ce qu'il a réellement obtenu. Le meilleur vaccin est celui avec le taux d'immunité le plus élevé, c'est-à-dire μ∗. Or le médecin peut avoir choisi des vaccins peut-être moins efficaces au cours des essais.
Note importante : pour la suite, les résultats seront généralement réalisés avec 100 initialisations différentes du MAB (tous les MAB ont 5 vaccins mais des taux d'immunistation différent) pour réduire le bruit de simulation. Concrètement, on exécutera au moins 100x generate_arms.
I. Cas classique des bandits manchots¶
I.a. Solution Gloutonne¶
Le médecin fonctionne sur deux phases :
- Exploration : Le médecin calcule le taux d'immunisation empirique sur les N premiers patients en administrant le même nombre de fois chaque vaccin :
$$\hat{\mu_i}[0\rightarrow N] = \frac{1}{T_i} \sum_{k=0}^{N-1} \chi_{v_k,i}R_k,$$
avec $T_i = \sum_{k=0}^{N-1} \chi_{v_k,i}$.
- Exploitation : Le vaccin $v_i = \arg\max_j \hat{\mu_j}[0\rightarrow N]$ est utilisé pour les M patients suivants.
Q2. Implémentez cette solution avec N = 50 et M = 500 et testez-la avec 100 MAB. On souhaite savoir si vous trouvez le vaccin optimal à l'issue d'une phase d'exploration. Quelle est l'espérance empirique de cette variable ? Et son écart-type ? Calculez de même l'espérance et l'écart-type du regret sur vos 100 simulations.
Pour rappel, le regret est défini par :
$$r_n = n\mu^* - \sum_{k=0}^{n-1} R_k$$
Attention : $n$ est le nombre total de patients, donc ici $N + M$.
N = 50 # Number of patients for exploration
M = 500 # Number of patients for exploitation
n = N + M # Total number of patients
def greedy_solution(MAB: t.List[ArmBernoulli], N: int = 50 , M: int = 500 ):
rewards = torch.zeros(K)
trials = torch.zeros(K)
# Exploration
for i in range(0,N,1):
arm_index = torch.randint(0, K, (1,)).item()
reward = MAB[arm_index].sample().float().item()
rewards[arm_index] += reward
trials[arm_index] += 1
empirical_means = rewards / trials
best_vaccine = torch.argmax(empirical_means).item()
# Exploitation
total_reward = rewards.sum().item()
for i in range(0, M, 1):
reward = MAB[best_vaccine].sample().float().item()
total_reward += reward
return total_reward, rewards, trials
# Simulation
def simulate(N: int = 50, M: int = 500, num_simulations: int = 100):
regrets = []
exploration_successes = []
for i in range(0, num_simulations, 1):
MAB = generate_arms(K) # K arms, K = 5
mu_star = find_mu_star_value(MAB)
total_reward, rewards, trials = greedy_solution(MAB, N, M)
total_patients = N + M
regret = total_patients * mu_star - total_reward
regrets.append(regret)
best_vaccine_index = torch.argmax(torch.tensor([arm.immunity_rate for arm in MAB])).item()
if (best_vaccine_index == torch.argmax(rewards / trials).item()):
exploration_successes.append(1)
else:
exploration_successes.append(0)
# Stats
mean_regret = np.mean(regrets)
std_regret = np.std(regrets)
mean_exploration_success = np.mean(exploration_successes)
std_exploration_success = np.std(exploration_successes)
return mean_exploration_success, std_exploration_success, mean_regret, std_regret
num_simulations = 100
mean_exploration_success, std_exploration_success, mean_regret, std_regret = simulate(N, M, num_simulations)
print(f"Espérance empirique de la réussite de l'exploration : {mean_exploration_success}")
print(f"Écart-type de la réussite de l'exploration : {std_exploration_success}")
print(f"Espérance empirique du regret : {mean_regret}")
print(f"Écart-type du regret : {std_regret}")
Espérance empirique de la réussite de l'exploration : 0.76 Écart-type de la réussite de l'exploration : 0.4270831300812525 Espérance empirique du regret : 33.71997833251953 Écart-type du regret : 40.201202392578125
Q3. On étudie maintenant l'influence de la taille du training set $N$. On considère que N+M est une constante, puis on fait varier N entre K et M. Calculez le regret pour ces différentes tailles du training set différents MAB et representez le regret moyen, le regret min et max (vous devriez trouver une courbe en U ou en V pour le regret moyen). Quelle est la taille optimale du training set ?
# Fonction pour simuler le regret pour une valeur donnée de N
def simulate_regret(N, M, MAB, num_simulations=100):
K = len(MAB)
total_patients = N + M
regrets_sim = []
for _ in range(num_simulations):
rewards = torch.zeros(K)
trials = torch.zeros(K)
total_reward = 0
# Phase d'exploration
for i in range(N):
arm_index = torch.randint(0, K, (1,)).item()
reward = MAB[arm_index].sample().float().item()
rewards[arm_index] += reward
trials[arm_index] += 1
total_reward += reward
best_vaccine = torch.argmax(rewards / trials).item()
# Phase d'exploitation
for i in range(M):
reward = MAB[best_vaccine].sample().float().item()
rewards[best_vaccine] += reward
total_reward += reward
regret = total_patients * find_mu_star_value(MAB) - total_reward
regrets_sim.append(regret)
# Même chose que pour simulate mais ici on stock les regrets
mean_regret = torch.tensor(regrets_sim).mean().item()
min_regret = torch.tensor(regrets_sim).min().item()
max_regret = torch.tensor(regrets_sim).max().item()
return mean_regret, min_regret, max_regret
def simulation_plot(num_simulations=100):
K = 5
N_plus_M = 550 # N + M constante
N_values = range(K, N_plus_M - K + 1, 50)
MAB = generate_arms(K)
mean_regrets = []
min_regrets = []
max_regrets = []
minimum_x = 0 # Coordonnées du minimum de regret
minimum_y = 0
minimum_value = float('inf')
# Pour chaque valeur de N
for N in N_values:
M = N_plus_M - N
mean_regret, min_regret, max_regret = simulate_regret(N, M, MAB, num_simulations)
# Get the minimum regret coordinates
if mean_regret < minimum_value:
minimum_value = mean_regret
minimum_x = N
minimum_y = mean_regret
mean_regrets.append(mean_regret)
min_regrets.append(min_regret)
max_regrets.append(max_regret)
return N_values, mean_regrets, min_regrets, max_regrets, minimum_x, minimum_y
# Visualiser les résultats
fig, axes = plt.subplots(1, 3, figsize=(18, 6))
simulations = [10, 50, 100]
for i, num_simulations in enumerate(simulations):
ax = axes[i]
N_values, mean_regrets, min_regrets, max_regrets, minimum_x, minimum_y = simulation_plot(num_simulations=num_simulations)
ax.plot(N_values, mean_regrets, label='Mean Regret', color='blue')
ax.fill_between(N_values, min_regrets, max_regrets, color='blue', alpha=0.2, label='Min/Max Regret')
ax.set_title(f"num_simulations={num_simulations}")
ax.set_xlabel("Taille du training set (N)")
ax.set_ylabel("Regret")
ax.legend()
ax.grid(True)
ax.plot(minimum_x, minimum_y, marker='o', color='red', markersize=5)
print(f'Le minimum de regret est atteint pour N={minimum_x} avec un regret de {minimum_y} concernant {num_simulations} simulations')
plt.tight_layout()
plt.show()
Le minimum de regret est atteint pour N=55 avec un regret de 25.3706111907959 concernant 10 simulations Le minimum de regret est atteint pour N=55 avec un regret de 13.861010551452637 concernant 50 simulations Le minimum de regret est atteint pour N=55 avec un regret de 35.238014221191406 concernant 100 simulations
On a bien la courbe en V pour le regret moyen. La taille optimale du training set est 55 pour 100 simulations. Au dela de 55, le regret moyen augmente. Cela veut dire que, passé ce seuil, le fait de continuuer à explorer appoorte plus de regret que de bénéfice.
Q4. On propose d'améliorer l'algorithme précédant en mettant à jour les taux d'immunisation empiriques $\hat{\mu}_i$ pendant la phase d'exploitation (algorithme greedy). Concrètement, à chaque nouvel patient, on lui administre le meilleur vaccin selon les stats. Notez vous une amélioration du regret ? Proposez un exemple où les taux d'immunisation du MAB ne changent rien.
def greedy_solution_with_update(N: int, M: int, MAB: t.List[ArmBernoulli]):
rewards = torch.zeros(K)
trials = torch.zeros(K)
# Exploration
for _ in range(N):
arm_index = torch.randint(0, K, (1,)).item()
reward = MAB[arm_index].sample().float().item()
rewards[arm_index] += reward
trials[arm_index] += 1
empirical_means = rewards / trials
best_vaccine = torch.argmax(empirical_means).item()
# Exploitation
total_reward = rewards.sum().item()
for _ in range(M):
reward = MAB[best_vaccine].sample().float().item()
total_reward += reward
rewards[best_vaccine] += reward
trials[best_vaccine] += 1
empirical_means[best_vaccine] = rewards[best_vaccine] / trials[best_vaccine]
return total_reward, rewards, trials
def simulate_with_update(N: int, M: int, num_simulations: int = 100):
regrets = []
exploration_successes = []
for _ in range(num_simulations):
MAB = generate_arms(K)
mu_star = find_mu_star_value(MAB)
total_reward, rewards, trials = greedy_solution_with_update(N, M, MAB)
# Regret
total_patients = N + M
regret = total_patients * mu_star - total_reward
regrets.append(regret)
# Vérification si le meilleur vaccin a été trouvé
best_vaccine_index = torch.argmax(torch.tensor([arm.immunity_rate for arm in MAB])).item()
if (best_vaccine_index == torch.argmax(rewards / trials).item()):
exploration_successes.append(1)
else:
exploration_successes.append(0)
# Stats
mean_regret = np.mean(regrets)
std_regret = np.std(regrets)
mean_exploration_success = np.mean(exploration_successes)
std_exploration_success = np.std(exploration_successes)
return mean_exploration_success, std_exploration_success, mean_regret, std_regret
N = 50
M = 500
mean_exploration_success, std_exploration_success, mean_regret, std_regret = simulate_with_update(N, M)
print(f"Espérance empirique de la réussite de l'exploration : {mean_exploration_success}")
print(f"Écart-type de la réussite de l'exploration : {std_exploration_success}")
print(f"Espérance empirique du regret : {mean_regret}")
print(f"Écart-type du regret : {std_regret}")
Espérance empirique de la réussite de l'exploration : 0.79 Écart-type de la réussite de l'exploration : 0.4073082370883261 Espérance empirique du regret : 34.90438461303711 Écart-type du regret : 49.04499816894531
On note une légère amélioration du reget. Un exemple où cette solution ne marcherait pas est : Si tous les vaccins ont des taux d'immunisation très proches les uns des autres, la mise à jour pendant l'exploitation pourrait ne pas changer le choix du meilleur vaccin. Dans ce cas, le regret resterait donc similaire à celui du glouton classique.
Q5. Nouvelle amélioration : à chaque nouveau patient, on choisit si on lui administre le meilleur vaccin avec une probabilité $\epsilon$ ou un vaccin aléatoire ($p=1-\epsilon$). Vérifiez si vous obtenez un meilleur résultat avec N = 0 ou N > 0. À votre avis, à quoi sert $\epsilon$ ?
def epsilon_greedy_solution(N, M, MAB, epsilon):
total_reward = 0
rewards = torch.zeros(len(MAB))
trials = torch.zeros(len(MAB))
# Phase d'exploration
for _ in range(N):
arm = random.randint(0, len(MAB) - 1)
reward = MAB[arm].sample().float().item()
total_reward += reward
rewards[arm] += reward
trials[arm] += 1
# Phase d'exploitation avec epsilon-greedy
for _ in range(M):
if random.random() < epsilon:
# Exploration: choisir un bras aléatoire
arm = random.randint(0, len(MAB) - 1)
else:
# Exploitation: choisir le meilleur bras connu
arm = torch.argmax(rewards / (trials + 1e-6)).item()
reward = MAB[arm].sample().float().item()
total_reward += reward
rewards[arm] += reward
trials[arm] += 1
return total_reward, rewards, trials
def simulate_epsilon_greedy(N, M, epsilon, num_simulations=1000):
regrets = []
exploration_successes = []
for _ in range(num_simulations):
MAB = generate_arms(K)
mu_star = max(arm.immunity_rate for arm in MAB)
total_reward, rewards, trials = epsilon_greedy_solution(N, M, MAB, epsilon)
total_patients = N + M
regret = total_patients * mu_star - total_reward
regrets.append(regret)
# Vérification si le meilleur vaccin a été trouvé
best_vaccine_index = torch.argmax(torch.tensor([arm.immunity_rate for arm in MAB])).item()
if (best_vaccine_index == torch.argmax(rewards / trials).item()):
exploration_successes.append(1)
else:
exploration_successes.append(0)
# Stats
mean_regret = np.mean(regrets)
std_regret = np.std(regrets)
mean_exploration_success = np.mean(exploration_successes)
std_exploration_success = np.std(exploration_successes)
return mean_exploration_success, std_exploration_success, mean_regret, std_regret
# Paramètres
N_values = [0, 10 , 25, 50]
M = 500
epsilon = 0.1
for N in N_values:
print(f"\nRésultats pour N = {N}:")
mean_exploration_success, std_exploration_success, mean_regret, std_regret = simulate_epsilon_greedy(N, M, epsilon)
print(f"Espérance empirique de la réussite de l'exploration : {mean_exploration_success}")
print(f"Écart-type de la réussite de l'exploration : {std_exploration_success}")
print(f"Espérance empirique du regret : {mean_regret}")
print(f"Écart-type du regret : {std_regret}")
Résultats pour N = 0: Espérance empirique de la réussite de l'exploration : 0.878 Écart-type de la réussite de l'exploration : 0.3272858078194042 Espérance empirique du regret : 34.29494857788086 Écart-type du regret : 24.850692749023438 Résultats pour N = 10: Espérance empirique de la réussite de l'exploration : 0.863 Écart-type de la réussite de l'exploration : 0.3438473498516456 Espérance empirique du regret : 28.681297302246094 Écart-type du regret : 18.347148895263672 Résultats pour N = 25: Espérance empirique de la réussite de l'exploration : 0.874 Écart-type de la réussite de l'exploration : 0.331849363416596 Espérance empirique du regret : 29.642946243286133 Écart-type du regret : 15.815997123718262 Résultats pour N = 50: Espérance empirique de la réussite de l'exploration : 0.89 Écart-type de la réussite de l'exploration : 0.31288975694324034 Espérance empirique du regret : 38.33185958862305 Écart-type du regret : 16.31777572631836
ma réponse, brute :¶
Epsilon sert à introduire le fait de pouvoir casser la routine et permettre à l'algorithme de ne pas s'engouffrer dans un seul chemin lors de son exploration. Imaginons qu'on ait plusieurs machines à sous et qu'on trouve dans nos premières expériences une qui fonctionne bien comparée aux autres. Alors dans notre entrain on ne va jouer que sur la même machine indéfiniment. Or on passe sûrement à côté d'autres occasion en ne gardant que la même option à chaque tour. Epsilon permet donc de temps en temps de pouvoir changer d'alternative et ainsi espérer obtenir de meilleurs gains.
ma réponse, reformulée par un LLM:¶
Rôle d'epsilon dans l'algorithme epsilon-greedy.
Le paramètre epsilon joue un rôle crucial dans l'algorithme epsilon-greedy en introduisant un équilibre entre l'exploitation et l'exploration. Voici comment il fonctionne :
Exploitation : La plupart du temps (1 - epsilon), l'algorithme choisit l'option qui semble la meilleure selon les informations actuelles.
Exploration : Une fraction epsilon du temps, l'algorithme choisit une option au hasard, indépendamment des performances passées.
Cette approche permet à l'algorithme d'éviter de se focaliser prématurément sur une seule option qui pourrait sembler optimale au début, mais qui pourrait ne pas l'être à long terme. Sans cette exploration, l'algorithme risquerait de passer à côté d'autres options potentiellement meilleures.
Prenons l'exemple des machines à sous :
- Sans epsilon : Si une machine donne de bons résultats au début, l'algorithme pourrait s'y tenir indéfiniment, ignorant potentiellement des machines plus rentables.
- Avec epsilon : Même après avoir trouvé une bonne machine, l'algorithme continuera occasionnellement à essayer d'autres options, augmentant ainsi ses chances de découvrir la meilleure machine à long terme.
En résumé, epsilon permet à l'algorithme de maintenir un équilibre entre l'exploitation des connaissances acquises et l'exploration de nouvelles possibilités, optimisant ainsi les gains potentiels sur le long terme.
I.b. Borne inférieure de Lai & Robbins [Lai et Robbins, 1985]¶
Lai et Robbins [Lai et Robbins, 1985] considère une classe d'algorithmes $\pi$ pour résoudre ce type de problèmes.
Ils ont trouvé une borne inférieure sur les récompenses cumulées en valeur asymptotique :
$$\lim_{n\rightarrow \infty} \inf_{\pi} \frac{\sum_{k=0}^{n-1} R_k}{\log n} \geq \sum_{i~\text{tel que}~\mu_i \lt \mu^*} \frac{\mu^∗−\mu_i}{\text{KL}(\mu_i, \mu^*)} :=C(\mu)$$
avec $\text{KL}(x, y) = x \log(x/y) + (1 − x) \log((1 − x)/(1 − y))$ (distance de Kullback-Leibler) et $\sum_{k=0}^{n-1} R_k$ la récompense obtenue sur $n$ patients.
Q6. Justifiez pourquoi on peut en déduire que le regret d'un algorithme raisonnable sera au pire logarithmique.
On cherche à minimiser le regret cumulé. Or, l'exploration diminue exponentiellement rapidement au fur et à mesure que 𝑛 n augmente, car le regret est proportionnel à $\log 𝑛$. En fait, à long terme, les bras sous-optimaux sont tirés de moins en moins souvent.
L'algorithme va trouver rapidement quel bras est optimal, c'est-à-dire le bras avec la plus grande $\mu_i$, et va commencer à l'exploiter. Les autres bras seront tirés de moins en moins souvent (et donc avec une fréquence exponentiellement décroissante) pour minimiser le regret.
Q7. Tracez le regret issu de la borne de Lai & Robbins et comparez le au regret obtenu avec l'algorithme glouton.
# Calcule de la distance KL
def kl_divergence(mu_i, mu_star):
return mu_i * np.log(mu_i / mu_star) + (1 - mu_i) * np.log((1 - mu_i) / (1 - mu_star))
# Calcule de la borne de Lai & Robbins
def lai_robbins_bound(mus, mu_star, n):
C_mu = sum((mu_star - mu_i) / kl_divergence(mu_i, mu_star) for mu_i in mus if mu_i < mu_star)
return C_mu * np.log(n)
def greedy_algorithm(MAB, N, M):
rewards = torch.zeros(K, dtype=torch.float)
trials = torch.zeros(K, dtype=torch.float)
total_reward = 0
regrets = []
n = N + M
# Exploration
for _ in range(N):
i = torch.randint(0, K, (1,)).item()
reward = MAB[i].sample().float().item()
rewards[i] += reward
trials[i] += 1
empirical_means = rewards / trials
best_vaccine = torch.argmax(empirical_means).item()
# Exploitation
total_reward = rewards.sum().item()
for _ in range(M):
reward = MAB[best_vaccine].sample().float().item() #
rewards[best_vaccine] += reward
trials[best_vaccine] += 1
total_reward += reward #
regret = n * find_mu_star_value(MAB) - total_reward
regrets.append(regret)
return regrets
K = 5
N = 50
M = 500
MAB = generate_arms(K)
mus = [arm.immunity_rate for arm in MAB]
mu_star = find_mu_star_value(MAB)
lai_robbins_regret = [lai_robbins_bound(mus, mu_star, n) for n in range(1, M + 1)]
greedy_regret = greedy_algorithm(MAB, N, M)
plt.figure(figsize=(12, 6))
plt.plot(lai_robbins_regret, label="Lai & Robbins Bound", color='red')
plt.plot(greedy_regret, label="Greedy Algorithm Regret", color='blue')
plt.xlabel("Number of patients (n)")
plt.ylabel("Cumulative Regret")
plt.title("Comparison of Lai & Robbins Bound and Greedy Algorithm Regret")
plt.legend()
plt.grid(True)
plt.show()
La borne inférieure de Lai et Robbins (1985) pour le problème des bandits multi-bras est donnée par :
$$ \lim \inf_{n \rightarrow \infty} \frac{\sum_{k=0}^{n-1} R_k}{\log n} \geq \sum_{i~\text{tel que}~\mu_i < \mu^*} \frac{\mu^* - \mu_i}{\text{KL}(\mu_i, \mu^*)} := C(\mu) $$
où
$$ \text{KL}(\mu_i, \mu^*) = \mu_i \log\left(\frac{\mu_i}{\mu^*}\right) + (1 - \mu_i) \log\left(\frac{1 - \mu_i}{1 - \mu^*}\right) $$
donc :
$$ \lim \inf_{n \rightarrow \infty} \frac{\sum_{k=0}^{n-1} R_k}{\log n} \geq C(\mu) $$
donc :
$$ \sum_{k=0}^{n-1} R_k \geq C(\mu) \cdot \log n $$
d'où :
$$ Regret_{Lai\&Robbins}(n) = C(\mu) \cdot \log(n) $$
Ici, on observe que le regret issu de la borne de Lai & Robbins est bien inférieur au regret obtenu avec l'algorithme glouton.
I.c. Upper Confidence Bounds¶
Cet algorithme améliore la version précédente en ajoutant un biais lié à la fréquentation de chaque vaccin :
$$\bar{\mu}_i = \hat{\mu}_i + \sqrt{\frac{C\log{n}}{T_i}}$$,
avec $C=2$.
Q8. Implémentez la modification de cet algorithme. Observez un intérêt à conserver $N > 0$ ? Et $\epsilon < 1$ ? Expliquez pourquoi.
Dans la suite, on prendra $N = 0$ et $\epsilon = 1$.
def ucb_algorithm(N, epsilon, C, num_vaccines, num_iterations):
successes = np.zeros(num_vaccines)
trials = np.zeros(num_vaccines)
total_reward = 0
rewards = []
for n in tqdm(range(1, num_iterations + 1)):
if n <= N:
vaccine = np.random.randint(num_vaccines)
else:
if np.random.random() < epsilon:
vaccine = np.argmax(successes / np.maximum(trials, 1) + np.sqrt((C * np.log(n)) / np.maximum(trials, 1)))
else:
vaccine = np.random.randint(num_vaccines)
reward = np.random.binomial(1, true_probs[vaccine])
total_reward += reward
rewards.append(total_reward / n)
successes[vaccine] += reward
trials[vaccine] += 1
return rewards, successes, trials
num_vaccines = 5
num_iterations = 100
true_probs = np.random.random(num_vaccines)
# Test 1: N = 0, epsilon = 1 (UCB pur)
rewards_ucb, successes_ucb, trials_ucb = ucb_algorithm(0, 1, 2, num_vaccines, num_iterations)
# Test 2: N > 0, epsilon = 1
rewards_N, successes_N, trials_N = ucb_algorithm(100, 1, 2, num_vaccines, num_iterations)
# Test 3: N = 0, epsilon < 1
rewards_eps, successes_eps, trials_eps = ucb_algorithm(0, 0.9, 2, num_vaccines, num_iterations)
plt.figure(figsize=(15, 5))
plt.plot(rewards_ucb, label='UCB pur (N=0, ε=1)')
plt.plot(rewards_N, label='N>0, ε=1')
plt.plot(rewards_eps, label='N=0, ε<1')
plt.xlabel('Itérations')
plt.ylabel('Récompense moyenne')
plt.legend()
plt.title('Comparaison des performances pour différentes configurations')
plt.show()
print("Taux d'immunisation empiriques:")
print("UCB pur:", successes_ucb / trials_ucb)
print("N>0:", successes_N / trials_N)
print("ε<1:", successes_eps / trials_eps)
print("\nVrais taux d'immunisation:", true_probs)
0%| | 0/100 [00:00<?, ?it/s]
0%| | 0/100 [00:00<?, ?it/s]
0%| | 0/100 [00:00<?, ?it/s]
Taux d'immunisation empiriques: UCB pur: [0. 0.55 0.3 0.82608696 0.55555556] N>0: [0. 0.57142857 0.28 0.78571429 0.75 ] ε<1: [0. 0.57142857 0.2 0.67567568 0.44444444] Vrais taux d'immunisation: [0.01738161 0.56080022 0.30341546 0.70008349 0.52924351]
- UCB pur (N=0, ε=1) : L'algorithme explore beaucoup (grâce au 𝜖=1), mais cela le ralentit au début, car il passe beaucoup de temps à tester des bras moins efficaces avant de converger vers les meilleurs.
- N>0, ε=1 : Cette configuration a un bon début, car elle fait une exploration initiale avec 𝑁>0 N>0, mais elle continue de beaucoup explorer, ce qui empêche une exploitation efficace des meilleurs bras dans les dernières itérations.
- N=0, ε<1 : Cet arrangement semble être le plus efficace, car il trouve rapidement les bons vaccins et les exploite bien. Ce qui conduit à une récompense moyenne plus élevée plus tôt dans les itérations.
Mais les résultats sont assez semblables donc tout cela reste à nuancer.
Q9. Tracez sous la forme d'une animation l'évolution des taux d'immunisation empirique (fig. de gauche) et l'évolution du regret (fig. droite). Dans la figure de gauche, vous representerez $\bar{\mu}_i$ et $\hat{\mu}_i$ pour chaque vaccin.
def ucb_algorithm(true_probs, num_iterations, C=2):
num_vaccines = len(true_probs)
empirical_means = np.zeros(num_vaccines)
counts = np.zeros(num_vaccines)
rewards = []
empirical_rates = []
regrets = []
for t in range(1, num_iterations + 1):
ucb_values = empirical_means + C * np.sqrt(np.log(t) / (counts + 1e-9))
chosen_vaccine = np.argmax(ucb_values)
# Simulate vaccination
reward = np.random.random() < true_probs[chosen_vaccine]
# Update statistics
counts[chosen_vaccine] += 1
empirical_means[chosen_vaccine] += (reward - empirical_means[chosen_vaccine]) / counts[chosen_vaccine]
rewards.append(reward)
empirical_rates.append(empirical_means.copy())
regrets.append(np.max(true_probs) - true_probs[chosen_vaccine])
return np.array(empirical_rates), np.cumsum(regrets)
# Set up the animation
num_iterations = 1000
#true_probs = np.array([0.3, 0.5, 0.7])
true_probs = np.random.random(5)
empirical_rates, cumulative_regrets = ucb_algorithm(true_probs, num_iterations)
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(15, 5))
# Lines for empirical rates (\hat{\mu}_i)
lines1 = [ax1.plot([], [], label=f'Vaccin {i+1}')[0] for i in range(len(true_probs))]
# Lines for theoretical true rates (\bar{\mu}_i)
true_lines1 = [ax1.axhline(y=prob, color=f'C{i}', linestyle='--', alpha=0.5) for i, prob in enumerate(true_probs)]
ax1.set_xlim(0, num_iterations)
ax1.set_ylim(0, 1)
ax1.set_xlabel('Itérations')
ax1.set_ylabel('Taux d\'immunisation')
ax1.set_title('Évolution des taux d\'immunisation empiriques')
ax1.legend()
# Line for cumulative regret
line2, = ax2.plot([], [], label='Regret cumulatif')
ax2.set_xlim(0, num_iterations)
ax2.set_ylim(0, np.max(cumulative_regrets) * 1.1)
ax2.set_xlabel('Itérations')
ax2.set_ylabel('Regret cumulatif')
ax2.set_title('Évolution du regret cumulatif')
ax2.legend()
def init():
for line in lines1:
line.set_data([], [])
line2.set_data([], [])
return lines1 + [line2]
def animate(frame):
for i, line in enumerate(lines1):
line.set_data(range(frame), empirical_rates[:frame, i])
line2.set_data(range(frame), cumulative_regrets[:frame])
return lines1 + [line2]
anim = FuncAnimation(fig, animate, frames=num_iterations, init_func=init, blit=True, interval=20)
HTML(anim.to_jshtml())
Animation size has reached 21052478 bytes, exceeding the limit of 20971520.0. If you're sure you want a larger animation embedded, set the animation.embed_limit rc parameter to a larger value (in MB). This and further frames will be dropped.